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We study the antiferromagnetic XY model on a triangular lattice by extensive Monte Carlo simulations, focusing on its 
ordering and critical properties. Our result clearly shows that two separate transitions occur at two distinct temperatures, 
the one at a higher temperature is associated with a Z2-symmetry breaking driven by the chirality, and the one at a 
lower temperature is associated with the onset of the quasi-long-range order of the XY spin. We carefully examine the 
critical properties of each transition to find that the criticality of the chiral transition is consistent with the standard two- 
dimensional Ising universality class, whereas that of the spin transition might differ from the conventional Kosterlitz- 
Thouless (KT) one. The observed non-KT nature of the spin criticality is consistent with the most recent simulation 
result on the fully-frustrated XY model on a square lattice. 
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1. Introduction 

The ordering of frustrated magnets has long been stud- 
ied quite extensively as one of primary subjects in statistical 
physics. '~ 5 ' Frustrated magnets show rich phase diagrams and 
new critical behaviors associated with novel physical objects 
which are composites of several spins. 3 ~ 6 ' A typical example 
might be the "chirality", which quantifies the handedness of 
the noncollinear spin structure. Recent theoretical and exper- 
imental studies have revealed that the chirality often plays an 
essential role in the ordering of many frustrated magnets. 5 ' 6 ' 

For last decades, various types of two-dimensional (2D) 
fully-frustrated XY (FFXY) models, which are believed to 
show a common critical behavior, have been intensively inves- 
tigated. 7 ~ 36 ' They include the 2D XY model on a square lattice 
with both the ferromagnetic and the antiferromagnetic (AF) 
couplings, 7 ~ 15 ' and that on a triangular lattice with the purely 
AF coupling. 7 ' 8 ' 13 ' 16 ^ 20 ' These 2D FFXY models are known 
to be characterized by the Z2 chiral degeneracy, i.e., the ex- 
istence of twofold degenerate ground states distinguished by 
the opposite signs of the chirality. The system then has a pos- 
sibility to exhibit a chirality-driven phase transition associated 
with the spontaneous breaking of this Z2 symmetry. 

Main motivation to study the 2D FFXY models has been to 
reveal whether the chiral transition occurs separately from the 
spin transition or not. After many controversial discussions 
about this point in past days, recent intensive studies have 
established that this "spin-chirality decoupling" phenomenon 
really occurs in the 2D FFXY models. Namely, two phase 
transitions take place at two distinct temperatures. 13 ~ 15 ' 21 ' The 
transition at a higher temperature T c involves the Z2 chiral- 
symmetry breaking, while the spin ordering at a lower tem- 
perature T s (< T c ) is accompanied by the onset of the quasi- 
long-range order of the spin similar to the Kosterlitz-Thouless 
(KT) transition. 37 ' 

Despite the progress concerning the existence of two sep- 
arate chirality and spin transitions, the critical properties of 
these transitions still remain somewhat controversial. For the 
chiral transition at T c , a plausible universality class might 
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be the 2D Ising one, since the chiral transition is accompa- 
nied by the Z2-symmetry breaking. Although several studies 
indeed supported the Ising universality, 1415 ' some counter- 
evidences, either based on theoretical considerations 36 ' or on 
numerical calculations, 13 ' 20 ' 26 ' 27 ' were also reported. 

Concerning the spin transition at T s < T c , the issue of 
whether it belongs to the standard KT universality class or 
not has also been hotly debated. Indeed, possible deviations 
of the critical exponents from the standard KT values were 
reported by several researches, 15 ' I6 ' 20 ' 26 ~ 28 ' 34 ' 35 ' while some 
others claimed that the spin transition belonged to the con- 
ventional KT universality class. 13 ' 14 ' 36 ' 

Under these circumstances, in order to obtain further in- 
sight into the critical properties of the chiral and the spin tran- 
sitions, we investigate in the present paper the ordering of the 
AF XY model on a triangular lattice by extensive Monte Carlo 
(MC) simulations. Among the two typical 2D FFXY mod- 
els, i.e., the square-lattice FFXY model and the triangular- 
lattice AF XY model, the former was recently studied quite 
extensively by large-scale equilibrium MC simulations up to 
the linear size as large as (9(1 3 ) by Hasenbusch, Pelissetto, 
and Vicari, 14 ' and by Okumura, Yoshino, and Kawamura. 15 ' 
Concerning the latter model, i.e., the triangular-lattice AF XY 
model, by contrast, a comparable numerical study has still 
been lacking. Equilibrium simulations were so far limited 
to the size of (9(1 2 ), e.g. up to L = 108 in ref. 20, while 
nonequilibrium simulations, which dealt with an apparently 
large size of L = 2000, had a limitation of the finite-time 
effect. 13 ' In order to fulfill this gap for the triangular-lattice 
AF XY model, and to get further information about its criti- 
cal properties, we perform here a large-scale equilibrium MC 
simulation on the model up to the size L = 1152, compara- 
ble to those recently performed for the square-lattice FFXY 
model. 

Our result suggests that the chiral transition belongs to the 
standard Ising universality, while the spin one belongs to the 
non-KT universality. This conclusion is fully consistent with 
that of the latest simulation on the FFXY model on the square 
lattice. 15 ' 

The remaining part of this paper is organized as follows. In 
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§2, we introduce our model and define several physical quan- 
tities, including the specific heat, the spin and the chiral order 
parameters (ordering susceptibilities), the spin and the chi- 
ral Binder ratios, the spin and the chiral correlation lengths, 
the helicity modulus and the vorticity modulus, which are 
computed in our simulations and employed in our analysis of 
the critical properties. Details of our MC simulations are also 
given here. In §3, we present the result of our MC simulations. 
The critical properties of both the chiral and the spin transi- 
tions are analyzed in detail. Finally, §4 is devoted to summary 
and discussion. 

2. MODEL AND PHYSICAL QUANTITIES 

Our model is the AF XY model on the triangular lattice. We 
assume that sites are labeled by i,j = 1,2, •• • ,N), and 
the corresponding coordinate is denoted as r, = (x,,y,)- The 
number of spins N is related to the linear system size L as 
N - Lx L. The XY spin on a site i, 5,, has two components 
as Si = (Sf,S>) = (cos U sin 0) where < < 2n. The 
Hamiltonian is given by 



(l) 



where J t j is the exchange interaction and kept constant J < 
for the rest of the paper, and the summation (ij) is taken over 
the all nearest neighbors. The partition function is given by 



-0F 



= Z = jf\%e-*", 
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where ft is the inverse temperature l/(k B T), £ B being the 
Boltzmann constant. In the following, we normalize the tem- 
perature in a unit of J/k B = 1. The thermal average will be 
denoted by the angular brackets (•••). 

2.1 Chirality '-related observables 

First, we define the local chirality on a plaquette p, an up- 
ward elementary triangle, as 
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where e z denotes a unit vector along the z direction. The di- 
rected sum Yj(i,j)€ P i s taken over three bonds surrounding a 
plaquette p in a counter-clockwise direction. The normaliza- 
tion factor is chosen to satisfy the condition \k p \ — 1 in the 
ground state where the spins take the so-called 120-degrees 
or the V3 X V3 structure. The chiralities exhibit a checker- 
board pattern in the ground state. The associated chiral order 
parameter is defined by, 
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where pj denotes the index of a plaquette whose bottom-left 
is the site j. 

The chiral Binder parameter is defined by 
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and the chiral correlation length is 
1 



(\m c \ 2 ) 



2sin(7r/L) V <|m c (qr min )| 2 ) 



-1, 
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where q m i n = (2n/L)e x . We introduce the Fourier component 
of the chirality with wavevector q by 



7=1 



The chiral susceptibility is defined by 

Xc = Np(mt). 



(7) 
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2.2 Spin- related observables 

The V3 x V3 spin structure can be captured via the TT-point 
magnetization of the wavevector q s = (4n/3, 0) 



1 N 

7=1 



The corresponding Binder parameter is defined by 

(Kl 4 ) 
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The spin correlation length is given by 
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where the Fourier transform of the y<f -point magnetization is 
given by 



1 N 

m s {q)=j i Yu e ' i(qs+q>riS i- 

7=1 
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The 7<f -point spin susceptibility is defined by 

Xs =Np(\m s \ 2 ). 

In addition to the above quantities, we calculate the helicity 
modulus to detect the rigidity of the system against a uniaxial 
twist applied to the spin. This quantity is known to effectively 
detect the KT transition. Consider an infinitesimal twist AG/L 
applied to the spin along, say, the x axis. In the Hamiltonian 
(1), we replace (0,- - 6j) by (0,- - 9j + (A6/L)e x ■ c,y), where dj 
is a lattice vector connecting the site i to the site j. Expanding 
the free energy with respect to A0, we get 

AF = Y'(A6) 2 + 6KA0 4 ). 



(14) 



where 



r = i{(Z- / ^-^ 2cos ^-^ 



-P 



^ Jij(e x ■ e u ) sin(0 ; - 6j) 

V<«,7> 



(15) 



Since the spin density of the triangular lattice is 2/ V3 times 
larger than the square lattice, 18 ' 20) we need to normalize Y' to 
make a direct comparison with the KT theory. The normalized 
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helicity modulus Y is defined by 

Y=^-Y'. (16) 
V3 

To obtain further insight into the low-temperature phase, 
we also calculate the vorticity modulus. 38 ' The vorticity mod- 
ulus measures the rigidity of the system against an emergence 
of an isolated vortex, and is originally calculated by the dif- 
ference between the total free energies under the two different 
boundary conditions (BC): One is the usual fixed BC where 
the boundary spins are fixed so that they match the ground- 
state V3x V3 spin configuration, and the other is the "vortex" 
BC where the boundary spins are fixed so that a single vortex 
is accommodated on top of the fixed-BC spin configuration. 38 ' 
More conveniently, the same quantity can be estimated by cal- 
culating appropriate fluctuations under the standard periodic 
BC, 15 ' 19) which we adopt in our present computation. Con- 
sider the phase shifts of spins caused by an emergence of an 
isolated vortex with a winding number m. The phase shifts are 
written as 



(Gj - Gj) -> (Gi - Gj + m<p u ) 



(17) 



where 



is an angle between the sites i and /' with 



respect to the vortex core. We assume that the vortex core is 
located at the center of the lattice. Conventionally, the solid 
angle is approximated by 0y = e-,j ■ ti/r] where r\ is the 
distance between site i and the vortex core, and ti is a unit 
vector enclosing the vortex core and tangential to the circular 
path passing through the site i. Neglecting the discreteness of 
the vortex winding number (m = 1 , 2, ■ • • ) and expanding the 
free energy with respect to m, we get the vorticity modulus V 
as 



where 



AF = Vm 2 + 0(m 4 ), 

V = lY Sitcom -0} 
■pi J^Jij^jSmiGi-Gj) 

U /; ;\ 



(18) 
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For large enough systems, the vorticity modulus is expected 
to behave as 



V(T,L) « c(T) + v(T)lnL, 



(20) 



where we call v(T) the reduced vorticity modulus. 

We expect that the present system gets the rigidity against 
a single vortex excitation at low temperatures, but loses it at 
high temperatures, i.e., v(T) is expected to take finite values 
at and below the spin transition temperature. This property 
might be utilized in estimating the spin transition temperature. 
In practice, we estimate v(T) by taking the difference between 
the vorticity modulus of different sizes as 

V(T,Li)-V(T,L 2 ) 



v(T\L u L 2 ) 



(21) 



ln(Li/L 2 ) 

for large L\ and L 2 . As an additional remark, we comment on 
the problem associated with the BC employed in computing 
the vorticity modulus. Under the periodic BC employed in the 
present simulation, <pjj becomes ill-defined at the boundary. 
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Fig. 1. (Color online) The temperature dependence of the specific heat per 
spin around the chiral transition temperature T c . The arrow indicates our es- 
timate of the transition temperature T c . 



We avoid this problem by setting 4>ij to be zero in computing 
the vorticity modulus if the bond between i and j passes across 
the boundary. 

2.3 Details of simulations 

The standard MC method with the Metropolis up- 
date combined with the over-relaxation update is used in 
our simulation. Our unit MC step per spin (MCS) con- 
sists of one Metropolis sweep followed by three con- 
secutive over-relaxation sweeps. System sizes are L = 
36,72, 144,288,576, 1152, and L = 48,96,192,384,768. 
Typical number of total MCSs are 4 x 10 6 MCS for L < 288, 
4 x 10 7 MCS for 384 <L< 576, and 10 8 MCS for 768 < L< 
1152. The initial 3 x 10 6 MCS are discarded for equilibra- 
tion. The BC is taken to be periodic. We extrapolate the data 
at the measuring temperature to nearby temperatures by us- 
ing the histogram method. This enables us to calculate phys- 
ical quantities in a wide temperature range by performing a 
small number of simulations, which is particularly useful in 
precisely locating the critical temperature. 

3. MONTE CARLO RESULTS 

In this section, we present the result of our MC simulations, 
and examine the critical properties of the chiral and the spin 
transitions. 

3.1 Chiral transition 

We first present in Fig. 1 the temperature dependence of 
the specific heat per spin, C(T) = f3 2 ((<K 2 ) - <-K> 2 ) /N. The 
data show a clear peak around the chiral transition tempera- 
ture T = T c . The size dependence of the peak temperatures 
Tpeak(L) is shown in Fig. 2. We see from the figure that there 
is a characteristic size around L x ~ 300 at which the size de- 
pendence of the peak location changes drastically. For smaller 
sizes L < L x , a near-linear dependence is observed, while, for 
larger sizes L> L x , another linear dependence with a different 
slope is observed. Similar crossover behavior was observed 
also for the square-lattice FFX Y model. 15) This crossover may 
be understood by considering the relation between the chiral 
and the spin correlation lengths around T c . 

In Fig. 3, we plot the finite-size spin correlation length % S {L) 
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Fig. 2. (Color online) The peak temperatures of the specific heat plotted 
versus the inverse linear size 1 /L. The upper figure shows the data for a wide 
range of L, while the lower one is a magnified view of larger sizes L > L x ss 
300. The dotted straight line (blue) is a linear fit of the data for L < 300, and 
the solid straight line (red) is the one for L > 300. The dotted curve (black) 
represents a nonlinear fit of the data of L > 300 on the basis of eq. (22) with 
v c a fitting parameter. 



versus the system size L at a temperature T = 0.51251 which 
is close to our estimate of T c (see below). One sees from the 
figure that the bulk spin-correlation length at the chiral transi- 
tion temperature is about £ s w 100, which is comparable to the 
observed crossover size L x * 300. Then, a natural explanation 
of the crossover behavior observed at L — L x might be the 
"spin-chirality decoupling" phenomenon occurring at longer 
length scales. Namely, the chiral correlation length outgrows 
the spin correlation length in the vicinity of T c at around 
the length scale L x . At shorter length scales or for smaller 
sizes, the spin and the chirality are trivially coupled, and the 
spin correlation apparently dominates the ordering. At longer 
length scale or for larger sizes, the chiral correlation exceeds 
the spin one around T c , and eventually dominates the ordering 
of the system. The observed size crossover might be regarded 
as a manifestation of this crossover from the spin-chirality 
coupling behavior for smaller sizes to the spin-chirality de- 
coupling behavior for larger sizes. This means that the asymp- 
totic critical behavior is expected only for larger systems of 
L > L x ~ 300. 

In view of this, we extrapolate the peak temperatures of the 
specific heat by using the data of L > L x only, on the basis of 
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Fig. 3. (Color online) The size dependence of the finite-size spin correla- 
tion length at T = 0.51251 =s T c . The value converges to f s ss 100 for large 
enough L. 



the asymptotic form, 

7peak(i) = Tc + (const)!- 1 '* 



(22) 



where v c is the chiral correlation-length exponent. The near- 
linear dependence observed in the lower panel of Fig. 2 im- 
plies v c — 1, suggesting the standard 2D Ising universality for 
the chiral transition. The straight-line fit with v c = 1 yields the 
chiral transition temperature r c _i s = 0.51251(1), as indicated 
by the straight line (red) in the lower panel of Fig. 2. If, on the 
other hand, we do not fix v c and leave it as a free parameter in 
the fit of eq. (22), we get T c = 0.51254(3) and v c = 0.72(25), 
which is represented by the dotted curve (black) in the lower 
panel of Fig. 2. The quoted error bars are estimated by the 
standard ^-square method. In fact, the resulting value of v c is 
still compatible with the Ising value v c = 1 within the error 
margin. 

Next, we analyze the temperature and size dependence of 
the chiral Binder parameter g c and the chiral correlation- 
length ratio f c /L, which are shown in Fig. 4 and Fig. 5, respec- 
tively. Both quantities exhibit the crossing behavior around 
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Fig. 4. (Color online) The temperature dependence of the chiral Binder 
parameter g c . 
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Fig. 5. (Color online) The temperature dependence of the chiral 
correlation-length ratio f c /L. 



the chiral transition temperature as expected. The size de- 
pendence of the crossing temperatures between the sizes L 
and 2L, T X {L), are shown in Fig. 6. Interestingly, the size de- 
pendence of T x changes drastically around the crossover size 
L x =i 300 as already observed above for the specific-heat- 
peak temperatures. Namely, for L < L x , T x tends to decrease 
with increasing L, while, for L > L x , it tends to increase with 
L. Asymptotic behavior is then expected only for larger sizes 
L> L x . We fit T x of both the chiral Binder parameter and the 
chiral correlation-length ratio in the range of L > 384 to the 
standard finite-size-scaling form with the correction term, 



T X (L) = T c + (const)L 



(23) 



where w is the leading correction-to-scaling exponent. The 
exponent a> should be common for both the Binder parameter 
and the correlation-length ratio, while the constant term is not 
so. 

We try the combined fit of all the three quantities investi- 
gated, i.e., the specific-heat peak temperatures, the crossing 
temperatures of the chiral Binder parameter, and those of the 
the chiral correlation-length ratio, with T c , v c , and u> c as fitting 
parameters. The result is T c = 0.51254(3), v c = 0.79(23), and 
a) — 1.2(1 1), which is our best unbiased estimate of T c , v c and 
u) c . The associated fitting curve is given in the middle panel 
of Fig. 6 (the solid curve in red). In the lower panel, we plot 
the total ^-square values of this fit against the r c -value as- 
sumed in the fit. If we instead fix the values of T c and v e to 
the Ising values as determined above from the specific heat, 
r c _i s = 0.51251 and v c = 1, the fit yields a> = 2.2(6), which 
is well consistent with the Ising value u — 2. The associated 
fitting curve is also given in the middle panel of Fig. 6 (the 
broken curve in magenta). 

For further information about the criticality of the chiral 
transition, we analyze the chiral susceptibility % c . In Fig. 7, 
we show the size dependence of Xc for several temperatures 
around T c , At the chiral transition temperature T c , we ex- 
pect a power-law behavior Xc 00 L 2 ~ nc for large enough sys- 
tems. A complete linearity is not observed at any temperature, 
which implies that the scaling region has not completely been 
reached. Hence, we analyze the data on the basis of the fol- 
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Fig. 6. (Color online) (Upper) The inverse size dependence of the crossing 
temperatures of the chiral Binder parameter and of the chiral correlation- 
length ratio. The crossing temperatures are the ones between sizes L and 2L. 
(Middle) A magnified view of the upper figure in the large-L region. The 
straight line (magenta) gives the fit of the specific-heat-peak temperatures 
with the Ising value v c = 1, which yields T c = 0.51251. The broken curves 
(magenta) are the combined fits of the chiral Binder parameter and of the chi- 
ral correlation-length derived under those T c and v c . The solid curves (red) 
are the unbiased combined fit of all three quantities without fixing any pa- 
rameters. (Lower) The total ^-square values of the combined fit of the three 
quantities plotted against the transition temperature T c assumed in the fit. The 
horizontal line represents the total ^-square value which is greater than the 
optimal value by unity, usually used as an error criterion. 



lowing scaling form with a constant correction term, 
Xc (T,L) = a(T) + b{T)L 2 -^ ( - T \ 



(24) 



By fitting the data for large sizes of L > 288, we obtain the 
effective value of the chiral anomalous-dimension exponent 
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Fig. 7. (Color online) (Upper) The chiral susceptibility divided by L 2 , 
Xc/L 2 , plotted versus the linear size L at several different temperatures around 
T c . The solid line (red) represents a fit to eq. (24) at T = 0.51251. The dotted 
line (black) expresses a X LT n with the Ising exponent Tj = 0.25, and is plotted 
as a guide to the eye. (Lower) The temperature dependence of the effective 
chiral anomalous-dimension exponent ^(T) evaluated through eq. (24). The 
estimated transition temperatures r c _i s = 0.51251(1) and T c = 0.51254(3) 
(see the main text) are given with the error bars. The Ising value rj = 1/4 is 
obtained around T = T c -i s = 0.51251. 



rfiF(T), and the results are plotted versus the assumed critical 
temperature in Fig. 7. We find that the Ising value tj — 1/4 is 
obtained around T = 0.51251, which is quite consistent with 
the T f -value deduced from the specific-heat-peak tempera- 
tures r c _i s = 0.51251(1) with assuming the Ising exponent 
v c = 1. Meanwhile, our unbiased estimate T c = 0.51254(3) 
associated with v c = 0.79(23) yields somewhat larger value 
for 77c, iff = 0.295(63), which is still compatible with the 2D 
Ising universality class within the error margin, though the 
error bar is now pretty large. 

Summarizing the above results, we find that our unbiased 
estimates of T c , v c and u> are compatible with the Ising uni- 
versality, though a new universality class is also admissible 
within, say, 20% accuracy for the exponents. If we presume 
the Ising value v c = 1, we can precisely locate the transition 
temperature T c = 0.51251(1). This value of T c shows a good 
agreement with the revalue derived from the chiral suscepti- 
bility under the Ising assumption r\ = 1/4. Such consistency 
certainly favors the Ising universality for the chiral transition, 
though further investigations are desirable to completely set- 
tle the issue. 



3.2 Spin transition 

Next, we move to the spin ordering at T — T s . The spin 
Binder parameter g s (T, L) and the spin-correlation-length ra- 
tio £ S (T,L)/L of the sizes L = 36,72,144,288,576,1152 
are shown in Fig. 8 and Fig. 9, respectively. Other sizes 
L = 48,96,192,384,768 are also computed, but we omit 
them here for clarity since they show a similar behavior. In 
Fig. 10, we also plot the reduced vorticity modulus v(T\2L, L) 
versus the temperature for L = 36, 72, 144, 288, 576. 
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Fig. 8. (Color online) The temperature dependence of the spin Binder pa- 
rameter g s for system sizes L = 36, 72, 144, 288, 576, 1 152. 




0.5035 0.504 0.5045 0.505 0.5055 0.506 0.5065 0.507 0.5075 

T 

Fig. 9. (Color online) The temperature dependence of the spin correlation- 
length ratio &/Z, for system sizes L = 36, 72, 144, 288, 576, 1 152. 



From Figs. 8-10, we can see immediately that the cross- 
ing temperatures for the spin related quantities are located 
in the temperature range significantly lower than the chiral 
transition temperature estimated above, T c 0.5125, thereby 
providing further evidence that the spin transition occurs at 
a temperature lower than the chiral transition, i.e., the the 
spin-chirality decoupling. We also see from these figures that 
the data at lower temperatures tend to exhibit a weaker and 
weaker splay-out behavior for larger sizes, and instead show 
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Fig. 10. (Color online) The temperature dependence of reduced vorticity 
modulus v(T\2L, L) with L = 36, 72, 144, 288, 576. 



a merging behavior. Such a merging behavior, typically ob- 
served in systems showing the KT transition, is certainly ex- 
pected for the present system exhibiting a spin quasi-long- 
range order. 

The size dependence of the crossing temperatures between 
the sizes L and 2L are summarized in Fig. 11. To estimate 
the spin transition T s , we fit the L-dependence of the crossing 
temperatures T X {L) to the form, 



T X (L) = T s + (const)(lnL)~ 



(25) 



The value of T s and the exponent x should be common for all 
the three quantities, while the constant factor is not so. 

We perform the two types of fits on the basis of eq. (25): 
In one, we fix x = 2, while, in the other, we treat x as a free 
fitting parameter. The former corresponds to the conventional 
KT fit, while the latter allows for the non-KT fit. The former 
KT fit with x = 2 yields TVkt = 0.5024(5), while the latter 
fit yields 7V nonKT = 0.5046(10) and x = 3.6(12). To see the 
relevancy of these two types of fits, we plot the associated x 2 - 
values of the fits as a function of the assumed revalue in the 
lower panel of Fig. 11, where the optimal x has been adjusted 
as a function of T s . As can be seen from the figure, the best fit 
is obtained for the non-KT fit with 7V non KT = 0.5046, but the 
KT-fit with 7Vkt = 0.5024 lies within the error margin due 
to the asymmetry of the x 2 -curve, and cannot be ruled out. 

In order to discriminate between the KT and the non- 
KT criticalities for the spin transition, we proceed to exam- 
ine other quantities. We first evaluate the spin anomalous- 
dimension exponent rj s by analyzing the /T-point spin suscep- 
tibility Xs- At and below the spin transition temperature T s , 
a power-law behavior Xs(T,L) oc L 2 -is(- r > is expected. In Fig. 
12, we plot the spin susceptibility Xs divided by L 2 as a func- 
tion of L on a log-log plot, and we can see a good linearity 
of the data, in which the scaling region is expected. Thus, the 
effective spin anomalous-dimension exponent rj e f(T) can be 
defined through the form 

X s (T,L) = a(T)L 2 -^ T \ (26) 

The estimated rff(T) is plotted in the lower panel of Fig. 
12. From the figure, we find rj s = 0.226(14) at 7V non KT = 
0.5046. If we assume the KT value of rj s = 0.25, we find 




0.005 0.01 0.015 0.02 0.025 0.03 0.035 0.04 0.045 

^2 



1/(ln Lf 




0.5 0.501 0.502 0.503 0.504 0.505 0.506 



Fig. 11. (Color online) (Upper) The crossing temperatures of the spin 
Binder parameter g s , of the spin correlation-length ratio /j s /L, and of the re- 
duced vorticity modulus v are plotted versus (lnL)~ 2 . The straight (magenta) 
and curved (red) lines are the fits on the basis of eq. (25). The straight ones 
correspond to the KT fit (x = 2 in eq. (25)) and the curves are the non-KT 
fit with treating x as a fitting parameter. T s and x are taken to be common 
in each fit. (Lower) The^ 2 values of the non-KT fit on the basis of (25) are 
plotted against the assumed Tj-values, where x is adjusted as a function of 
T s . The horizontal straight line denotes the minimum ^ 2 -value plus unity. 



T s = 0.5063(1), which happens to be significantly higher than 
?Vkt = 0.5024(5) estimated above from the crossing temper- 
atures of the Binder ratio, of the correlation-length ratio and 
of the vorticity modulus. This discrepancy indicates that the 
KT criticality rj s = 0.25 and x = 2 might be inappropriate in 
describing the spin transition of the present model. 

We also employ the helicity modulus to further exam- 
ine the consistency. In 2D XY models, the helicity modu- 
lus is expected to show a discontinuous jump of magnitude 
T s I (Inrj S {T s )) at the transition temperature T s , This quantity 
relates the transition temperature T s to the exponent tj s . The 
temperature dependence of the helicity modulus is shown in 
Fig. 13. In the figure, we show the two straight lines T/(2m] s ), 
one with the KT value of rj s - 1/4 and the other with the 
non-KT value rj s = 0.226 as estimated above. By extrapo- 
lating the crossing temperatures between these two straight 
lines and the helicity-modulus data to L — oo on the basis 
of (25), we can obtain an estimate of T s . The resulting ex- 
trapolation curves are also given in Fig. 13. For the KT case 
77, = 1/4 and x = 2, we find T s = 0.5051(1), which, however, 
turns out to be inconsistent again with the revalues estimated 
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Fig. 12. (Color online) (Upper) The system-size L dependence of the spin 
susceptibility divided by L 2 at several temperatures on a log-log plot. The 
solid lines demonstrate the power law ZT''* with Tj s = 0.25 and r] s = 0.226. 
The conventional KT value r\ s = 0.25 is denoted by the broken line (green) 
obtained by fitting the data at T = 0.5063, the temperature associated with 
the KT transition, while the solid line (red) represents rj s = 0.226 obtained 
by fitting the data at T = 0.5046, the temperature associated with the non-KT 
transition (see the main text). Note that the Xs(T) data of L = 768 and 1152 
in a temperature range T < 0.5038 is obtained by a polynomial extrapola- 
tion of Xs(T) at higher temperatures. (Lower) The effective spin anomalous- 
dimension exponent rj s plotted versus the assumed T,-value. The estimated 
transition temperature Tj-nonKT = 0.5046(10) is also given with the error 
bar. The KT transition temperature Tj-kt = 0.5024(5) is located outside this 
temperature range. 



above by other independent methods under the KT assump- 
tion, i.e., T s = 0.5024(5) and T s = 0.5063(1). On the other 
hand, for the non-KT case rj s = 0.226 and x = 3.6, we get 
T s = 0.5044(1) which turns out to accord well with the T s - 
value, 7V non KT = 0.5046(10), estimated above by the other 
methods. We also try the fit where x is treated as a free fit- 
ting parameter. The result is T s = 0.5043(1) and x = 3.2(18), 
which yields essentially the same T s as the constrained non- 
KT fit. 

Considering all these results, we may now exclude the stan- 
dard KT criticality for the spin transition. The spin transition 
of the present model has a non-KT character. Our best un- 
biased estimates of the spin-transition temperature and the 
spin anomalous dimension are T s = 0.5046(10) and rj s - 
0.226(14), respectively. We note that the possibility of the 
non-KT nature of the spin transition was suggested before 
by several authors. 16 ' 20 ' Here we demonstrate it with a higher 



Fig. 13. (Color online) (Upper) The temperature dependence of the helicity 
modulus of different sizes. The lower solid straight line (red) and the upper 
broken straight line (green) denote T/(2m] s ) with = 1/4 and r) s = 0.226, 
respectively. (Lower) The intersection points of the helicity-modulus data and 
the straight lines of the upper figure are plotted against (In L)~ 2 . The solid line 
(red) is the fit under the KT assumption T] s = 1/4, x = 2, and the broken line 
(green) is the fit under the non-KT assumption = 0.226, x = 3.6. The 
former KT fit yields T s = 0.5051(1), whereas the latter non-KT fit yields 
T s = 0.5044(1). The broken line (blue) is the fit with T] s = 0.226 and with x 
as a free fitting parameter. 



confidence level by performing larger-scale simulations. We 
should also mention that the recent result on the square-lattice 
FFXY model by Okumura et a/ 15) gave essentially the same 
conclusion. The observed consistency between the square and 
triangular lattices may suggest a generality of the non-KT na- 
ture of the spin transition in the 2D FFXY models. 

4. Summary and Discussions 

In this paper, we investigated the AF XY model on the tri- 
angular lattice by extensive MC simulations, focusing on its 
spin and chirality orderings. Large systems up to L = 1152 
were simulated. Our data unambiguously showed that the chi- 
ral and spin transitions took place at two different tempera- 
tures, i.e. the occurrence of the spin-chirality decoupling. The 
spin transition temperature turned out to be lower than the 
chiral transition temperature by about 1.5%. A crossover be- 
havior from the spin-chirality coupling regime to the decou- 
pling regime was observed around the size L x m 300. This 
crossover is associated with the interchange of the relative 
dominance of the spin and the chiral correlation lengths. 

The universality class of the chiral and the spin transitions 
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was also examined. Concerning the criticality of the chiral 
transition, although we cannot completely rule out the pos- 
sibility of the non-Ising universality class, the high level of 
consistency we observed among different independent ob- 
servables certainly speaks for the Ising universality class. If 
we presume the Ising nature of the chiral transition, i.e., pre- 
sume v c = 1 , the chiral transition temperature was estimated 
quite accurately as T c = 0.51251(1). If we do not assume the 
Ising nature of the chiral transition, we get T c = 0.51254(3), 
v c = 0.79(23) and rj c = 0.295(63) which is still compatible 
with the Ising criticality, though the error bars are pretty large. 

Let us now compare our present estimates of T c , v c , and 
j] c with the recent estimates by other authors for the same 
model. 13 ' 20 ' On the basis of an equilibrium MC simulation, 
S. Lee and K.-C. Lee gave T c = 0.513(1), v e = 0.833(8) 
and T] c = 0.25(3), 20) which accord with ours. However, their 
data were taken only for rather small sizes of up to L = 108, 
smaller than the crossover length scale L x ~ 300, so that 
the data might fail to capture the genuine scaling behavior. 
As such, the accordance of their results with our present esti- 
mates may partly be accidental. On the basis of the nonequi- 
librium relaxation method, Ozeki and Ito gave T c = 0.512(1), 
v c = 0.84(2) and r\ c = 0.25(1), claiming a new universal- 
ity class for the chiral transition. 13) The method employed, 
i.e., the nonequilibrium relaxation method, however, has an 
intrinsic limitation of the finite observation time in spite of 
the apparently large system size studied. 

For the FFXY model on the square lattice, a model different 
from the present one but believed to lie in a common univer- 
sality class, two extensive equilibrium MC simulations com- 
parable in their size to our present one, have recently been 
made. Hasenbusch et al. studied the sizes up to L = 1000, 
and concluded that the chiral transition belonged to the 2D 
Ising universality class. 14 ' More precisely, these authors esti- 
mated the effective critical exponents for finite-size systems 
and observed a tendency converging to the Ising criticality for 
larger sizes of L > 500. Interestingly, a qualitative change of 
the size dependence was observed around the crossover size 
L x ~ 500, which accords with our present observation for the 
triangular XY model. Okumura et al. performed a large-scale 
MC simulation of the same model up to the size L = 1024. 15) 
These authors also concluded the Ising universality for the 
chiral transition, corroborating the conclusion of Hasenbusch 
et al. Thus, there now seems to be a reasonable amount of 
numerical evidence that the chiral transition of the 2D FFXY 
model lies in the standard 2D Ising universality class. 

Concerning the criticality of the spin transition, we ob- 
served an unacceptable level of inconsistency between several 
independent quantities studied when we assumed the stan- 
dard KT universality class for the spin transition. Such an in- 
consistency was resolved once we assumed the non-KT crit- 
ical exponents for the spin transition. Then, we concluded 
that the spin transition of the present model belonged to the 
non-KT universality class. As our final estimates, we quote 
T s = 0.5046(10) and t] s = 0.226(14). 

Comparing our present estimates of T s and T] s with the 
recent estimates by other authors on the same model, we 
find that our values accord with neither T s = 0.501(1) and 
r] s = 0.22(1) by Lee et al., 39) nor T s = 0.508(1) by Ozeki et 
a/. 13) Concerning the former case, the deviation of T s is pre- 
sumably due to the smallness of the simulation sizes in ref. 20. 



In view of the discrepancy in the estimated revalues, how- 
ever, the agreement in the estimated // s -values might proba- 
bly be accidental. Concerning the latter work by Ozeki et al, 
the cause of the discrepancy might come from the tendency 
that the short-time observation regime usually employed in 
the nonequilibrium relaxation method overestimates the tran- 
sition temperature. 

We also compare our estimate of rj s with the recent esti- 
mates by Hasenbusch et al. 14) and by Okumura et al. ly> for 
the square-lattice FFXY model. Hasenbusch et al. claimed 
the KT nature of the spin transition, based on their obser- 
vation that the two independent quantities, i.e., the helicity 
modulus and the spin correlation-length ratio, gave consistent 
results if one assumed the KT nature of the transition. Here 
we wish to comment that the type of their analysis was only 
a consistency check, which in itself is not enough to exclude 
the possible non-KT character of the transition. Okumura et 
al. examined the possibility of the non-KT nature of the spin 
transition on the same model, and reported the ^-value as 
T] s = 0.20(1), which was slightly smaller than our present es- 
timate T] s = 0.226(14), though they come near the quoted error 
margin. The difference of rj s between the two FFXY models, 
if it really occurs, would imply the non-universal character of 
the criticality of the spin transition of the 2D FFXY models. 
Further study is required to fully substantiate such a possibil- 
ity, however. 

The Ising universality of the chiral transition suggested in 
this paper is seemingly quite natural, since the chirality is 
an Ising variable possessing a Z2 symmetry. Nevertheless, 
we stress that the criticality of the chiral transition is still a 
nontrivial issue. In fact, in some systems such as vector spin 
glasses, critical exponents of the chiral transition are known 
to be quite different from the corresponding Ising values (i.e., 
from those of the Ising spin glass). 6 ' 40 ' 41) Hence, a naive sym- 
metry consideration is not always enough in identifying the 
universality class. A possible cause of such a deviation from 
the symmetry-based prediction might be the long-range na- 
ture of the chirality-chirality interaction. In 2D XY models, 
the effective interaction between the chiralities is known to be 
Coulombic, i.e., logarithmic. In the present 2D FFXY model, 
such a long-range nature of the chirality interaction appar- 
ently does not affect the universality class. This may be due 
to the simplicity of the relevant phase space of the present 
model, in contrast to that of spin glasses. It will then be inter- 
esting to clarify the conditions of the chiral transition belong- 
ing to the standard Ising universality class. 

The universality class of the spin transition is also a quite 
nontrivial issue. Because of the spin-chirality decoupling, the 
spin transition of the present model occurs in the presence of 
the chiral order already established at a higher temperature 
T c . This implies that the phase space available for the spin 
degrees of freedom at T s < T c is rather restricted. Such a 
situation is very different from the spin transition of the con- 
ventional unfrsutrated 2D XY models where an entire phase 
space is available, and hence, the nature of the spin order- 
ing of the current system could deviate from the conventional 
KT transition. To understand the universality of the transition, 
several theoretical approaches were already proposed, such 
as the Coulomb-gas representation of FFXY models 15 ' and 
the renormalization-group analyses of the coupled XY mod- 
els which are believed to be in the same universality class as 
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the 2D FFXY models. 34 ~ 36) However, a full comprehension of 
the phenomenon is still lacking. Another approach based on 
the KT renormalization group was proposed, 28 ) which yielded 
the nonuniversal jump of the helicity modulus, although the 
connection to the 2D FFXY models is not necessarily clear. 
Thus, further studies are required to fully understand the na- 
ture of the criticality of the spin transition of the 2D FFXY 
models. 
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for Scientific Research on Priority Areas 'Novel State of Mat- 
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1) G. H. Wannier: Phys. Rev. 79 (1950) 357. 

2) I. Shoji: Prog. Theor. Phys. 6 (1951) 306. 

3) P. W. Anderson: Mater. Res. Bull. 8 (1973) 153. 

4) J. Villain: J. Phys. C 10 (1977) 4793. 

5) H. Kawamura: J. Phys. Condens. Matter 10 (1998) 4707. 

6) H. Kawamura: J. Phys. Soc. Jpn. 79 (2010) 01 1007. 

7) J. M. Thijssen and H. J. F. Knops: Phys. Rev. B 42 (1990) 2438. 

8) J. Lee, J. M. Kosterlitz, and E. Granato: Phys. Rev. B 43 (1991) 11531. 

9) S. Lee and K.-C. Lee: Phys. Rev. B 49 (1994) 15184. 

10) G. Ramirez-Santiago and J. V. Jose: Phys. Rev. Lett. 68 (1992) 1224; 
Phys. Rev. B 49 (1994) 9567. 

1 1) P. Olsson: Phys. Rev. Lett. 75 2758 (1995); Phys. Rev. B 55 (1997) 3585. 

12) E. H. Boubcheur and H. T. Diep: Phys. Rev. B 58 (1998) 5163. 

13) Y. Ozeki and N. Ito: Phys. Rev. B 68 (2003) 054414. 

14) M. Hasenbusch, A. Pelissetto, and E. Vicari: Phys. Rev. B 72 (2005) 
184502; M. Hasenbusch, A. Pelissetto, and E. Vicari: J. Stat. Mech. 



(2005) P12002. 

15) S. Okumura, H. Yoshino, and H. Kawamura: Phys. Rev. B 83 (2011) 
094429. 

16) S. Teitel and C. Jayaprakash: Phys. Rev. B 27 (1983) 598. 

17) S. Miyashita and H. Shiba: J. Phys. Soc. Jpn. 53 (1984) 1 145. 

18) D. H. Lee, J. D. Joannopoulos, J. W. Negele, and D. P. Landau: Phys. 
Rev. Lett. 52 (1984) 433; Phys. Rev. B 33 (1986) 450. 

19) H.-J. Xu and B. W. Southern: J. Phys. A 29 (1996) L133. 

20) S. Lee and K.-C. Lee: Phys. Rev. B 57 (1998) 8472. 

21) S. E. Korshunov: Phys. Rev. Lett. 88 (2002) 167007. 

22) J. Lee, E. Granato, and J. M. Kosterlitz: Phys. Rev. B 44 (1991) 4819. 

23) E. Granato and M. P. Nightingale: Phys. Rev. B 48 (1993) 7438. 

24) M. P. Nightingale, E. Granato, and J. M. Kosterlitz: Phys. Rev. B 52 
(1995)7402. 

25) J. M. Thijssen and H. J. F. Knops: Phys. Rev. B 37 (1988) 7738. 

26) G. S. Grest: Phys. Rev. B 39 (1989) 9267. 

27) J.-R. Lee: Phys. Rev. B 49 (1994) 3317. 

28) P. Minnhagen: Phys. Rev. Lett. 54, 2351 (1985); P. Minnhagen: Phys. 
Rev. B 32, 7548 (1985). 

29) P. Minnhagen and M. Wallin: Phys. Rev. B 36 (1987) 5670; P. 
Minnhagen and M. Wallin: Phys. Rev. B 40 (1989) 5109. 

30) E. Granato and D. Dominguez: Phys. Rev. B 71 (2005) 094521. 

31) P. Olsson and S. Teitel: Phys. Rev. B 71 (2005) 104423. 

32) L. Capriotti, R. Vaia, A. Cuccoli, and V. Tognetti: Phys. Rev. B 58 (1998) 
273. 

33) S. Lee, K.-C. Lee, and J. M. Kosterlitz: Phys. Rev. B 56 (1997) 340. 

34) M. Y. Choi and D. Stroud: Phys. Rev. B 32 (1985) 5773. 

35) E. Granato, J. M. Kosterlitz, and J. Poulter: Phys. Rev. B 33 (1986) 4767. 

36) G. S. Jeon, S. Y. Park, and M. Y. Choi: Phys. Rev. B 55 (1997) 14088. 

37) J. M. Kosterlitz and D. J. Thouless: J. Phys. C 6 (1973) 1181. 

38) H. Kawamura and M. Kikuchi: Phys. Rev. B 47 (1993) 1 134. 

39) We estimate this value of rj s from the jump size of the helicity modulus 
given in ref. 20. 

40) H. Kawamura and M. S. Li: Phys. Rev. Lett. 29 (2001) 187204. 

41) D. X. Viet and H. Kawamura: Phys. Rev. Lett. 102(2009) 027202; D. X. 
Viet and H. Kawamura: Phys. Rev. B 80 (2009) 064418. 



10 



